% /*******************************************************************************
%  * Main Simulation Routine
%  * *****************************************************************************
%  * @param lower
%  * @param upper
%  * @param m
%  * @param n
%  * @param r
%  * @param k
%  * @param q
%  * @param t
%  * @param z
%  */
function [max_err_s, max_err_o, max_err_v, max_err_h] = get_max_err_all(lower,upper,m,n,r,k)

step_a = 2^(-m);
step_b = 2^(-n);

max_err_s = 0.0;
max_err_o = 0.0;
max_err_v = 0.0;
max_err_h = 0.0;


step = 2^(-(n+k));
f_t = 0:step:step*50;
e_star_f = zeros(size(f_t));

for i=1:length(f_t)
    e_star_f(i) = get_max_err(lower, upper, m, n, r, k, f_t(i));
end

max_err_s = min(e_star_f);


for a = lower:step_a:upper
    for b = lower:step_b:upper

        ans_true = a * b;
                
        [ans_const, corr_const] = ccm(a, b, m, n, r, k);
        abs_err_o = abs(ans_const - ans_true);
        if (abs_err_o > max_err_o)
            max_err_o = abs_err_o;
        end
        
        [ans_eric, corr_eric] = vcm(a, b, m, n, r, k);
        
        abs_err_v = abs(ans_eric - ans_true);
        if (abs_err_v > max_err_v)
            max_err_v = abs_err_v;
        end
        
        
        [ans_hyb, corr_hyb] = hcm(a, b, m, n, r, k, 0.5);
        
        abs_err_h = abs(ans_hyb - ans_true);
        if(abs_err_h > max_err_h)
            max_err_h = abs_err_h;
        end
            
    end
end
end